#include <math.h>
#include "fobjects.h"
#include "remez.h"

static int CalcN(const BandInfo &a, const BandInfo &p)
{
  Float loga = log10(a.GetErrorRatio()), logp = log10(p.GetErrorRatio());
  Float B = fabs(a.GetHigh() - p.GetLow()) / (2 * M_PI);
  Float F = 0.51244 * (logp - loga) + 11.012;
  Float D = ((0.005309 * logp + 0.07114) * logp - 0.4761) * loga - ((0.00266 * logp + 0.5941) * logp - 0.4278);
  return ((int)(sqr(D - F * B) / B + 1.5));
}

int DesignMultiBandFIR(const Filter F, Float d, FilterTypes ft, SearchTypes st, Vector &H)
{
  int i, j, N;
  Vector V;
  for(N = 0, i = 1; i < GetBandsNo(F); i++)
  {
    j = CalcN(F[i], F[i + 1]);
    if(j > N)
      N = j;
  }
  if(ft == OddAntysymetrical || ft == OddSymetrical)
  {
    if(!(N % 2))
      N++;
  }
  else
    if(N % 2)
      N++;
  InitRemez(F, N, ft);
  if(Remez(st, 20) != 0)
    return(-1);
  if(fabs(Delta) > d)
    for(N += 2; fabs(Delta) > d; N += 2)
    {
      InitRemez(F, N, ft);
      if(Remez(st, 20) != 0)
        return(-1);
      V = CalcH();
    }
  else
    for(N -= 2; fabs(Delta) <= d; N -= 2)
    {
      V = CalcH();
      InitRemez(F, N, ft);
      if(Remez(st, 20) != 0)
        return(-1);
    }
  H = V;
  return (0);
}

int DesignDifferentiator(Float UpLimit, Float ER, Float d, SearchTypes st, Vector &H)
{
  const int N1 = 5, N2 = 7;
  int N;
  Filter F(1);
  Float Delta1, Delta2;
  Vector V;
  F[1] = BandInfo(0, UpLimit, ER, 1);
  InitRemez(F, N1, OddDifferentiator);
  if(Remez(st, 20) != 0)
    return (-1);
  Delta1 = fabs(Delta);
  V = CalcH();
  InitRemez(F, N2, OddDifferentiator);
  if(Remez(st, 20) != 0)
    return (-1);
  Delta2 = fabs(Delta);
  if(d >= Delta2 && d < Delta1)
  {
    H = V;
    return (0);
  }
  N = (int)(N1 + log(d / Delta1) / log(d / Delta2) * (N2 - N1) + 0.5);
  if(N % 2 == 0)
    N++;
  InitRemez(F, N, OddDifferentiator);
  if(Remez(st, 20) != 0)
    return (-1);
  if(fabs(Delta) > d)
    for(N += 2; fabs(Delta) > d; N += 2)
    {
      InitRemez(F, N, OddDifferentiator);
      if(Remez(st, 20) != 0)
        return(-1);
      V = CalcH();
    }
  else
    for(N -= 2; fabs(Delta) <= d; N -= 2)
    {
      V = CalcH();
      InitRemez(F, N, OddDifferentiator);
      if(Remez(st, 20) != 0)
        return(-1);
    }
  H = V;
  return (0);
}
